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Abstract 

We study an autocatalytic system consisting of several interacting chemical species. We observe a 
strong dependence of the concentrations of the chemicals on the size of the system. This dependence is 
caused by the discrete nature of the molecular concentrations. Two basic mechanisms responsible for them 
are identified and elucidated. The relevance of the transitions to processes in biochemical systems and in 
micro-reactors is briefly discussed. 
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Rate equations are often employed in the study of biochemical reaction processes. In rate equations, 
the quantities of chemicals are treated as continuous variables, and the actual discreteness of the molecular 
concentration is ignored. Of course, fluctuations of numbers of molecules have been studied using stochastic 
differential equations and introduced non-trivial effects jl], |^|. Still, discreteness has not been considered in 
any such study. In many biochemical processes, however, some chemicals play important roles at extremely 
low concentrations, amounting to only a few molecules per cell 0]. Furthermore, there exist amplification 
mechanisms involving enzymes in cells through which even a change by one molecule in a cell can result in 
drastic effects. In such situations, the discreteness of the molecular concentration is obviously not negligible. 

We previously showed the existence of a novel transition induced by the discreteness of the molecular 
concentration in an autocatalytic reaction system Q. The system contains four chemicals Xi(i — 1, • • • ,4). 
We considered an autocatalytic reaction network (loop) represented by X- L + Xi + i — > 2Xi + i (with X§ = X±) 
within a container that is in contact with a reservoir of molecules. Through interaction with the reservoir, each 
molecule species Xi diffuses in and out at a total rate of DiSiV, where Di is the flow rate, Si the concentration 
of chemical Xi in the reservoir, and V the volume of the container. In this system, a novel state appears as 
a result of fluctuations and the discreteness of the molecular concentration, characterized as extinction and 
subsequent reemergence of molecule species alternately in the autocatalytic reaction loop. 

When the volume of the container is small, A^, the number of molecules of species Xj, may go to (i.e. 
become extinct) through a finite-size fluctuation due to the discreteness of the molecular concentration. Once 
Ni reaches 0, it remains until an Xi molecule flows in. Thus, if the flow rate of molecules is sufficiently small, 
state with N± = JV3 = or JV2 = N4, = can be realized. In a state with N2,N^ w (a "1-3 rich" state), 
switches between states with Ni > N3 and N3 > Ni can occur, and similarly for a state with N±,N3 w Q (a 
"2-4 rich" state). A symmetry-breaking transition to these states was observed in our previous study, with the 
decrease of V, as is shown in Fig. |], as the change of the probability distribution of the number of molecules. 
For large V, corresponding to the continuum limit, the distribution of z = ((N\ + N3) — (N2 + N4)) shows a 
single-peaked distribution around 2 = 0, whereas it is replaced by a symmetric, double-peak distribution as 
V is decreased. This is a novel discreteness-induced transition (DIT) occurring with the decrease of V. The 
transition occurs without any change of parameters, and thus cannot be discussed in the rate equation with 
noise (i.e., by the continuum description). 

In the system investigated in our previous work, the long-term average concentration of each chemical 
does not differ from that in the continuum limit, since the system over time switches between the two states 
of broken symmetry. It is important to determine if there are systems for which the average concentrations of 
chemicals are significantly altered by the DIT. We will show that this is possible in a system possessing some 
kind of asymmetry. With the DIT without symmetry breaking, the average concentrations of chemicals are 
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drastically altered by the change of V. The peak position of the distribution is changed with a finite jump, as 
the volume is decreased (see Fig. |2|, while see later sections for the description of the model and simulation). 
Borrowing the term of thermodynamics, the DIT reported previously is regarded as a second order transition 
involving symmetry breaking, while the DIT reported here corresponds to the first order transition without 
symmetry breaking. This result is biologically significant as providing a possible description of the alteration 
of the concentrations of some molecules within cells. Note in a cell, the number of molecules of each chemical 
species is not necessarily huge, and the discreteness effect is not always negligible. 

To investigate this problem, we again use an autocatalytic reaction loop of chemicals, but here we consider 
the case in which Di, r{ or Si is dependent on the chemical species i, where fj is the reaction constant of the 
reaction Xi + X i+ i — » 2X i+ i. To study the effects of discreteness, we investigated the reaction model by using 
a stochastic particle simulation. We assumed that the chemicals are well stirred in the container. At each 
simulation step, two molecules in the container are randomly chosen. Then we judge if the molecules react or 
not, by checking if one of the two acts as a catalyst for the other as a substrate. To carry out the simulation 
efficiently here, we adopted Gillespie's direct method || (see Appendix) Qj. 

Note that for there to be a DIT to a 1-3 or 2-4 rich state, it is necessary that the time interval for inflow of 
molecules be longer than the time scale of the reactions. This time interval for Xi inflow should be ~ 1/DiSiV. 
In our previous study, in which we considered the case of identical parameter values for all i, the discreteness 
of the molecular concentration has the same effect for all the molecule species, and the transition occurs near 
DsV = r. 

It is important to realize that the relevance of the discreteness of the molecular concentration depends 
on the reaction and flow rates of each molecule when the parameter values are not identical. For example, 
if DiSi < D 2 s 2 , the inflow time interval for X x molecules is longer than that for X 2 molecules, so that the 
discreteness of the Xi flow has a greater effect on the behavior of the system. In general, V may determine 
the species that become extinct, and the average concentration of each molecule can be greatly changed by 
the discreteness effect. 

Here, we consider the case in which Sj is species dependent, while D and r are identical for all species, 
for the autocatalytic loop introduced above. With this choice, the discreteness effect of each chemical Xi 
depends on i. In this system, we find discreteness effects that result in changes of the average concentrations 
Xi, with the temporal average of the concentration Xi = Ni/V. Although this result is obtained with this 
simple example, the mechanism we find would appear to be quite general, and hence there is reason to believe 
that the DIT we find exists in a wide variety of real systems. 

We first consider the effect of the discreteness of the inflow of chemicals and how this depends on the 
relation between the reaction rate and the inflow rate. In our model, the inflow interval of X% is ~ g 1 v , and 
the time scale of the reaction is ~ — — — . When the former time scale is larger than the latter, the reaction 
from i to i + 1 can proceed to completion before the inflow of species i occurs. Then Ni becomes 0. As long 
as Ni = 0, no reaction to produce chemical Xi occurs, and the average density may be decreased radically 
from the continuum limit case. 



Case I : inflow discreteness and reaction rate 

As a simplest example to study this mechanism, we consider the case with s\ = S3 > s 2 = S4. In this case, 
the rate equation in the continuum limit has a stable fixed point Vi : Xi — s^. When V is large, each x^ 
fluctuates around this fixed point. The average concentration x~i is shown in Fig. ^. With the decrease of V, 
the difference between the pair X\ and X3 and the pair x 2 and X4 is amplified, and there is clear deviation 
from the continuum limit case. 

The mechanism responsible for this amplification can be understood as follows. As V is decreased, we 
have found that the 1-3 rich state, with extinction of N 2 and A/4, appears when the increase of x~\ and X3 
occurs. To realize N 2 = N4 = 0, it is necessary for the inflow interval of X 2 or X4 to be longer than the time 
scale of the reaction. The inflow interval of Xi molecules is ~ g~y , while the time scale for the reaction is 
~ — — — . Since we set r = 1 and Xi = O(l), the 1-3 rich state appears for „ 1 , „ 1 > -, while the 2-4 
rich state appears for D ^ v , D l v > 7- Thus in the present case, the 1-3 rich state is first observed as V is 

decreased. In the range of values of V for which the relations D l v , D l 3 y < p < D l 2 y , D l iV are satisfied, 
the 1-3 rich state is realized often, while the 2-4 rich state is not (see Fig. |]). 

Once the 1-3 rich state is realized, an X 2 molecule and an X4 molecule must enter the system almost 
simultaneously for the system to break out of this state. Thus the 'rate of interruption' of the 1-3 rich state is 
roughly proportional to s 2 saV 2 , the product of the rates of X 2 inflow and X4 inflow. The expected residence 
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time in the 1-3 rich state is the reciprocal of the rate of interruption. Thus the ratio of the expected residence 
times in the 1-3 rich and 2-4 rich states is 2ila, 

S2S4 

From the above considerations, we expect that for some V « tt^i there appears a transition to the 1-3 
rich state, leading to a drastic increase of the 1-3 concentration. The validity of this conclusion has been 
confirmed by several simulations, one of whose results is shown in Fig. ^. 

Case I' : imbalance of inflow discreteness 

The transition discussed above can create a stronger effect on the concentrations. As an example, consider 
the case s\ = S3 > S2 > S4. 

In this case, as in case I, the 1-3 rich state is stable. While in this state, the system switches from a 
condition of N% > N3 to one of Ai < A3 due to A2 inflow and from Ni < A3 to Ai > A3 due to A4 inflow. 
Since the latter event is less frequent for S2 > S4, the condition iVi < A3 is satisfied for a greater amount 
of time in the 1-3 rich state. Hence, it is expected that x~i < £3. This is confirmed by the results displayed 
111 Fig. |. This is in strong contrast with the result in the continuum limit, where 11 k i 3 if D < rsi (i.e., 
the time scale of the reactions is much shorter than that of the inflow). The significant difference between x~\ 
and £3 found here appears only when the 1-3 rich state is realized through the effect of the discreteness of 
the flow of A4 molecules. As shown in Fig. ||, there is amplification of the difference between £\ and £3 as V 
decreases that occurs simultaneously with the transition to the 1-3 rich state. 

Case II : inflow and outflow 

When DsiV is small enough to insure the existence of both 1-3 and 2-4 rich states, the preference of states 
can depend on the concentrations s^. The preferred state is selected through another DIT caused by outflow 
rather than inflow of a particular chemical. 

As an example, we consider the case S2 > si > S3 = S4. Here again, the rate equation in the continuum 
limit has a stable fixed point. If D <C rsj, then x\,x% « Sl l" s? and X2,x i w 22±2± a t the fixed point. 

As discussed above, 1-3 and 2-4 rich states appear for small DS3V . In the 2-4 rich state, it is likely for A4 
to decrease as a result of the outflow of A4 and the reaction X4 + X\ — > 2 Ai facilitated by the inflow of X\ . 

If S4V < 1, it may be the case that all A4 molecules flow out, and N4 becomes 0. The time required to 
realize A4 = from N4 = n should be ~ logn when S4V w and n 3> 1. However, if s\ is large, X4 will 
be consumed by the reaction caused by X\, and for this reason, N4 will decrease to more rapidly. The time 
required to use up X4 may also depend on si. In this case, the 1-3 rich state is favoured by the mechanism 
described below. 

When A4 > 0, N4 can increase again as a result of A3 inflow, which leads to switch from the N2 > N4 
condition to the N2 < N4 condition. The inflow interval for A3 is ~ d1 3 v • ^ ^ ms m t erva l is much shorter 
than the time required for N4 to reach 0, N4 may increase again, causing the 2-4 rich state to be preserved. 
However, if the interval is longer, A4 may decrease to 0, in which case, the 2-4 rich state can be readily 
destroyed by the inflow of an X3 molecule, as shown in Fig. ^. 

When the system is in the 1-3 rich state, on the other hand, switches from a condition of N\ > N3 to one 
of Ni < N3 due to X2 inflow, cause A3 to remain large (as in case I'). The system therefore tends to maintain 
the condition Ni < A3 (as long as D/rsi is not too large). However, N± only rarely decreases to 0, unlike A4, 
because si is relatively large. Thus the 1-3 rich state is more stable than the 2-4 rich state, as shown in Figs. 
[7] and |. 

Hence, when V is decreased sufficiently to satisfy S4V < 1, and the time interval D v is sufficiently long 
to allow A4 to decrease to 0, the 2-4 rich state loses stability, and the residence time in the 1-3 rich state 
increases, due to the discreteness of A 4 . As a consequence, X3 increases as V decreases, as shown in Fig. ||. 

Amplification by Discreteness 

Summarizing the findings discussed above, differences among the 'degrees of discreteness' of the chemicals 
lead to novel DIT. The average chemical concentrations are greatly altered by this DIT. Indeed, as the system 
size (the volume V) changes, there is a sharp transition to a state qualitatively different from that found in 
the continuum limit. There are two key parameters with regard to discreteness: One is D l v (investigated in 
case I), the inflow time interval for Xi, and the other is siV (investigated in case II), the number of species 
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Xi molecules in the system when it is at equilibrium with the reservoir. 

If the interval n 1 is longer than the time scale of the reaction, — — — , the discreteness of the Xi inflow 
is relevant. In such a situation, the Xi molecules present in the system may be completely consumed by the 
reaction before any new Xi molecules flow in, so that N may reach 0. 

Then, if the condition SiV < 1 is satisfied in addition to the above stated condition, Ni can become as 
a result of all Xi molecules flowing out of the system. In this case, the relation between the time necessary 
to realize a switch that increases Ni and the time necessary for Ni to decay to is also important. 

With the above two conditions satisfied for each species Xi, there appear several switches to different 
states as V is changed. As an example, we considered the case in which si = 0.09, s 2 = 3.89, S3 = S4 = 0.01, 
and D = 1/64. In this case, the average concentration x~i exhibits three transitions as V is decreased, as 
shown in Fig. |To| . 

First, in the continuum limit, X\ and aT 3 are very small, as resulted from the fact that si + s 3 <C s 2 + s 4 . 
Around V = 10 3 , the discreteness of A3 becomes significant, and the 2-4 rich state appears. Then the reactions 
X-2 + A3 — > 2A3 and A3 + X4 — » 2A4 take place only sporadically. Contrastingly, the flow of X\ molecules 
is fairly steady. Thus, while the system is in the 2-4 rich state, N2 > A4 is satisfied for most of the time, as 
shown in case I'. Figure |ll| displays the distribution of X2- Double peaks corresponding to the 2-4 rich state 
appear in this situation. 

In the 2-4 rich state with N2 > N4,, A2 molecules flowing into the system raise N2 to the level of S2V, 
establishing equilibrium with the reservoir. At the same time, A4 molecules flow out, and A4 decreases to 
the level of S4V, as seen in case II. 

As seen in case II, The difference between N 2 and A4 increases with further decrease of V, since the 
switching rate decreases. In Fig. O, the gap between the two peaks in the distribution of A 2 is seen to 
become larger as V decreases. Around V — 10 2 , finally, the imbalance between N 2 and N 4 destabilizes the 
2-4 rich state. For this reason, the 1-3 rich state becomes almost as stable as (or more stable than) the 2-4 
rich state, in spite of the relation s% + S3 >C s 2 + S4. The residence time in the 1-3 rich state increases sharply, 
causing X3 to increase (as shown in Figs. ^ and |l0| ). In fact, 3T3 increases to approximately 2, which is more 
than 30 times larger than its value in the continuum limit || . 

For very small V (i.e. V < 2), N\ and A3 decrease to quite readily, and thus the 1-3 rich state is also 
easily destroyed. In this situation, for most of the time only one chemical exists in the container. Here, only 
x~2 has a large value, with all of the others near or at Ni = 0. 

In the manner described above, non-trivial alteration of chemical concentrations as a result of DIT was 
observed. It has been found that those molecule species whose numbers vanish are determined not only by 
the flow rates but also by the network and dynamics of the reactions. For example, when V is relatively large 
(V w 10 2 ), x~s decreases as S3 increases. 



Discussion 

In conclusion, we have reported a DIT that leads to a strong effect on the average concentrations of the 
chemicals. Although we have studied a simple case with only four chemicals here, we have found that this 
type of DIT appears in more complex reaction networks of a more general nature. 

In fact, we have randomly chosen a catalytic reaction network consisting of few hundred species, and 
studied the population dynamics of each chemical species with the scheme of stochastic simulation adopted 
here. For some reaction networks we have observed the DIT as the volume is decreased. In such cases, we 
have found that the combination of the two mechanisms we studied here leads to a variety of transitions and 
alterations of molecular concentrations. Although the example reported in the present paper is quite simple 
and may look special, the mechanism found in the example gives a basis for DIT in complex reaction network. 

Generally speaking, DIT and its effects on molecular concentrations are likely to be observed with chemical 
networks containing autocatalytic reactions. However, in some examples they are observed even without 
autocatalytic reactions. When some part of reaction networks works as an autocatalytic sub-network as a set, 
as seen in hypercycles |llj , the DIT of the present mechanism is possible. 

It is now experimentally feasible to construct a catalytic reaction system in a micro-reactor, and to design 
other types of systems with small numbers of molecules. Also, there was great advance in techniques for 
detection of small numbers (on the order of 1 to 10 2 ) of molecules using fluorescence or other new methods 
such as thermal-lens microscopy |Q. In such systems, experimental verification of DIT should be possible. 
Also, we believe that the alterations of chemical concentrations resulting from DIT that we found will have 
practical applications, since quite high accumulation of dilute chemical species is possible as we have shown 
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here. 

Since the number of molecules in a biological cell is often small, the relevance of DIT to cell biology is 
obvious. For example, in cell transduction, the number of signal molecules is often less than 100, and even a 
single molecule can switch the biochemical state of a cell Jll| . In our visual system, a single photon in retina is 
amplified to a macroscopic level 14|] . Transmission of signals through neurons via synapses also often involves 



a small number of molecules jl5|. Chemical reaction network consisting of several autocatalytic reaction is 
widely seen in a cell, and such autocatalytic process provides a candidate for amplification of an effect of a 
single molecule. Since the DIT we reported here is generally observed in autocatalytic reaction networks, it 
is expected that it may be used in a biochemical reaction network in a cell. Indeed, according to our results, 
the non-trivial accumulation of dilute molecules and switching among several distinct states with different 
chemical compositions may be realizable by, for example, the control of flow by receptor. Additionally, in 
some preliminary simulations with large reaction networks, some sub-networks can be effectively activated or 
inactivated by DIT. In such cases, transitions between several states characterized by active sub-networks can 
be observed, which will be relevant to switching between cellular states by a few signal molecules. 

Switching the expression of genes on and off is a focus of interest in bioinformatics. This digital behavior 
is also connected with the concentration of proteins present. As is pointed out Jl^| , genetic regulation is under 
stochasticity coming from smallness in the number of associated molecules. As we have seen in our model, 
one chemical species can exhibit both an on/off switch and continuous regulation of other chemicals, even if 
the number of molecules of this species is small. We believe that the switching of chemical states facilitated 
by our DIT plays a role in the regulation of genetic and metabolic processes in cells. 

Throughout the paper we have adopted stochastic particle simulations. Of course master equation ap- 
proach is also equivalently possible, which is especially useful if some analytic tools for it are developed. For 
example, use of Fokker-Planck equations derived in the limit of large volume (molecule numbers), is a powerful 
tool |L6| . Since our DIT occurs when the volume (the number of molecules) is quite small, such tools are so 
far not available. In future it will also be important to develop some analytic tools for a system where the 
discreteness in the number is essential. 
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Appendix : Details of the Simulation 

Details of the Model 

We assumed that the chemicals are well stirred in the container and the molecules have no volume. Thus 
the rate of the reaction X, + Xj_|_i — > 2Xi + i is given by Ri = riXiXi + \ [concentration / time], where is 
the reaction constant and Xi is the concentration of the chemical Xi. By rewriting it with the use of iVj, the 
number of Xi molecules, the rate of the reaction is given by RiV = r ' jv '^' +1 [reactions / time]. In the same 
way, the rate of the Xi inflow is given by DiSiV [molecules / time], corresponding to DiSi [concentration / 
time], whereas that of the Xi outflow is given by DiNi [molecules / time], corresponding to DiXi [concentration 
/ time]. 

This stochastic model approaches the rate equation 

when one takes a continuum limit, given by V — > oo. 

We also assumed that the area of the surface of the container is proportional to the volume V, and thus the 
rate of the Xi flow is proportional to DiV , to have this well-defined continuum limit for V — > oo. One might 
assume that the area of the surface should be V 2 ^ 3 , and the rate of the Xi flow should be proportional to 
DiV 2 / 3 . This change of setting alters just the parameter values. By suitably adjusting parameters Di and/or 
Sj, the same transitions to the switching states and the alteration of average concentrations are observed, even 
with such settings. 

With the rates of the reactions and the flows above, we carried out the stochastic simulation. In principle, 
one can carry out the simulation, by randomly selecting two molecules, and transforming one of them to other 
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molecule, according to the reaction rule, with the probability proportional to the rate of reaction, when these 
molecules react. Here, as an efficient simulation method, we adopt Gillespie's direct method, instead. 

Gillespie's Direct Method 

In our system, the state of the system is determined by iVj, the number of molecules, and is changed only 
when one reaction or one molecular flow occurs. Thus the rate of the reactions and the flows do not change 
until the next event (one reaction or one molecular flow) occurs, so that the lapse time to the next event 
decays exponentially. 

Gillespie's direct method |(J stands on this fact. First, we determine the lapse time to the next event by 
exponentially-distributed random numbers, and set the time forward. Next, we determine which event occurs, 
with the proportion to the rate of the event. We change the state according to the event, and re-calculate the 
rate of the reactions and the flows. These steps are executed repeatedly, until the specified time elapses. 

In some cases, especially with complicated reaction networks, there are more efficient methods (see refs. 
|l7| , Here, for simplicity, we adopted the Gillespie's direct method. Our result discussed above does not 
depend on which method to use. 
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Figure 1: Probability distribution of z = (x\ + X3) — (x 2 + X4). Here Si = 1, and D = 1/128. For V > 256, 
z has a distribution around 0, corresponding to the fixed point state Xi — 1 (for all i). For V < 64, the 
distribution has double peaks around z = 4, corresponding to the 1-3 rich state (Ni, N3 3> N 2l N4 w 0), and 
z = —4, corresponding to the 2-4 rich state. 




Figure 2: Probability distribution of X3 in our model with si = sq, = 1.99, S3 — S4 — 0.01, D = 1/128, 
sampled over a time span of 5 x 10 6 , for different V (see Case II for details). 




Figure 3: The average concentration Xi in Case I : Sx = S3 = 1.7, s 2 = S4 = 0.3, as a function of the volume 
V (sampled over a time span of 10 6 for V > 1024, 10 7 for 32 < V < 1024, and 10 8 for V < 32, also the same 
for Figs. I and |). D = 1/128. 
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Figure 4: Probability distribution of x\ and X2, sampled over a time span of 5 x 10 6 . si = S3 = 1.7, 
S2 = S4 = 0.3 (Case I). D = 1/128. In the case V = 2048, peaks around x, = Sj, which correspond to the 
fixed point at the continuum limit, is shown. As V decreases, the peaks get broader according to fluctuations, 
and the tail of the distribution of xi reaches 0. Thus, there appears a peak around X2 = 0. 
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Figure 6: Time series of the number of molecules Ni for V = 16. Here si = 0.09, S2 — 3.89, S3 = S4 = 0.01, 
and D = 1/64. There is a transition to the 2-4 rich state at t — 4511. In the 2-4 rich state, X4 molecules 
flow out at the rate D, and N4 thereby decreases. Due to the flow of X3 molecules, switching from N2 > N4 
to A^2 < N4 occurs, and N4 increases again (as seen at t = 4743). Here, the interval over which the switching 
takes place is longer than the interval of X3 inflow, l/DVs^ = 400, and is indeed long enough for most X4 
molecules to diffuse out before X3 molecules can flow into the system. Thus, here N4 readily decreases to 
before the switch. At t — 5300, with N4 — an X3 molecule flows into the container, leading to a switch to 
the 1-3 rich state. 
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Figure 7: Probability distribution of Xi in Case II : Si = s% = 1.99, S3 = S4 = 0.01, sampled over a time 
span of 5 x 10 6 , for different V. D = 1/128. In the case V — 512, all Xi shows peaks around Xi = 1, which 
correspond to the fixed point at the continuum limit. When V is small, each x\ shows a peak at a different 
concentration, and the peak height changes greatly with change of V. For V — 256, there appear peaks at 
around x\ = 1, X3 = 3, X2 = 4, and also around Xi = (for all i). With further decrease of V, N2 and 7V4 
reach more easily than N\ and N3 do. Accordingly the peaks for xi and X4 around grows as shown in the 
case V = 128. The system tends to stay at the 1-3 rich state, and the peak around x-i = 4 gets smaller (as 
shown in the case V = 64). 




Figure 8: Probability distribution of (x% + X3) and (x% + £4) in Case II : si = S2 = 1.99, S3 = S4 = 0.01, 
sampled over a time span of 5 x 10 6 . V = 32, D = 1/128. With such a small V, the 2-4 rich states are 
destabilized, and the rate of the residence at the 1-3 rich states is almost 10 2 times larger than that at the 
2-4 rich states. The state allowed by the continuum limit, X\ + X3 « X2 + X4 is very rare. 




Figure 9: The average concentration Xi in Case II : Si = S2 = 1.99, S3 = S4 = 0.01, as a function of the 
volume V. D = 1/128. 
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Figure 10: The average concentrations xt, sampled over a time span of 5 x 10 8 for V < 32 and 5 x 10 6 for 
V > 32, plotted as functions of the volume V. si = 0.09, s 2 = 3.89, s 3 = s 4 = 0.01, D = 1/64. For large V, 
x~i is close to the fixed point value of the continuum limit. As V decreases, there first appears a 2-4 rich state, 
but for smaller V, this 2-4 rich state becomes unstable, and the residence time of the 1-3 rich state increases, 
leading to a sharp increase of aT 3 . For much smaller V (< 0.5), only X 2 molecules exist for most of the time. 




Figure 11: Probability distribution of x 2 , sampled over a time span of 5 x 10 6 . s% = 0.09, s 2 = 3.89, 
S3 = S4 = 0.01, and D = 1/64. When V is large, there is a single peak around x 2 = 2, which corresponds to 
the fixed point in the continuum limit. Around V — 10 3 , double peaks appear around x 2 = 1 and x 2 = 3, 
corresponding to the 2-4 rich state. As V decreases, these two peaks move apart, and near V = 10 2 , the tail 
of the lower peak reaches 0, implying that N 2 and JV4 often decrease to 0. Thus here, the 2-4 rich state is 
unstable. Then, as V is decreased further, the peak at x 2 = rises sharply. 
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